(1) Load packages and read data

The dataset is downloaded from Metabolights. Accession number: MTBLS487. The data was obtained in positive ion mode over an m/z range of 300–2000 using a Thermo MALDI LTQ Orbitrap XL instrument (Thermo Fisher Scientific, Bremen, Germany). The spatial resolution was 50 micrometers. For more information about the dataset, readers are kindly referred to the original publication (Hall et al., 2017).

library(Cardinal)
library(BiocParallel)
library(plotly)
path <- paste0("Data/test_POS", ".imzML")
Brain <- readMSIData(path, resolution = 5, units = "ppm", 
                     mass.range = c(300, 2000), attach.only = T)

(2) Pre-processing

Brain2 <- 
  Brain %>%
  normalize(method = "tic") %>%
  peakPick(method = "mad", SNR = 6) %>%
  peakAlign(tolerance = 5, units = "ppm") %>%
  process(BPPARAM = SerialParam())

Check the ion image of m/z 844.527 before and after pre-processing in order to make sure that chosen pre-processing parameters are fine.

darkmode()
##(2.1) before pre-processing
image(Brain, mz = 844.527, smooth.image = "gaussian", plusminus = 0.03, colorscale = magma,
      contrast.enhance = "suppression", normalize.image = "linear")

##(2.2) after pre-processing
image(Brain2, mz = 844.527, smooth.image = "gaussian", plusminus = 0.03, colorscale = magma,
      contrast.enhance = "suppression", normalize.image = "linear")

(3) Co-localization

Here first 100 colocalized ions are retained.

#(3.1) co-localization PC36:1
coloc_826 <- colocalized(Brain2, mz= 826.572, n = 100, BPPARAM = SerialParam())
#(3.2) co-localization PC38:6
coloc_844 <- colocalized(Brain2, mz= 844.525, n = 100, BPPARAM = SerialParam())
#(3.3) co-localization PC40:6
coloc_872 <- colocalized(Brain2, mz= 872.557, n = 100, BPPARAM = SerialParam())

(4) Extract in-source fragments, and create pseudo MS/MS spectrum

#(4.1) PC36:1
mycoloc <- as.data.frame(coloc_826)
MSe <- mycoloc[mycoloc$correlation >= 0.90,]
## extract intensity
int = as.vector(rep(NA, dim(MSe)[1]))
for (i in 1:length(MSe$mz)) {
  int[i] <- sum(spectra(Brain2)[features(Brain2,  mz = MSe$mz[i]),])
}
## interactive plot
spec = cbind.data.frame(mz = MSe$mz, Int = int)
p = ggplot(spec, aes(x = mz, ymax = Int/max(Int)*100, ymin = 0, colour = "red")) +
  geom_linerange() +
  scale_y_continuous(expand = c(0, 0), limits = c(0, 100*1.1)) +
  theme_bw() +
  ggtitle("m/z 826.577")
ggplotly(p)
#(4.2) PC38:6
mycoloc <- as.data.frame(coloc_844)
MSe <- mycoloc[mycoloc$correlation >= 0.90,]
## extract intensity
int = as.vector(rep(NA, dim(MSe)[1]))
for (i in 1:length(MSe$mz)) {
  int[i] <- sum(spectra(Brain2)[features(Brain2,  mz = MSe$mz[i]),])
}
## interactive plot
spec = cbind.data.frame(mz = MSe$mz, Int = int)
p = ggplot(spec, aes(x = mz, ymax = Int/max(Int)*100, ymin = 0, colour = "red")) +
  geom_linerange() +
  scale_y_continuous(expand = c(0, 0), limits = c(0, 100*1.1)) +
  theme_bw() +
  ggtitle("m/z 844.528")
ggplotly(p)
#(4.3) PC40:6
mycoloc <- as.data.frame(coloc_872)
MSe <- mycoloc[mycoloc$correlation >= 0.90,]
## extract intensity
int = as.vector(rep(NA, dim(MSe)[1]))
for (i in 1:length(MSe$mz)) {
  int[i] <- sum(spectra(Brain2)[features(Brain2,  mz = MSe$mz[i]),])
}
## interactive plot
spec = cbind.data.frame(mz = MSe$mz, Int = int)
p = ggplot(spec, aes(x = mz, ymax = Int/max(Int)*100, ymin = 0, colour = "red")) +
  geom_linerange() +
  scale_y_continuous(expand = c(0, 0), limits = c(0, 100*1.1)) +
  theme_bw()+
  ggtitle("m/z 872.559")
ggplotly(p)

(5) Plot MALDI images

#(5.1) PC36:1
mycoloc <- as.data.frame(coloc_826)
MSe <- mycoloc[mycoloc$correlation >= 0.90,]
pdf(file = file.path("Result/colocalization_826.pdf"), onefile = TRUE)
for(i in 1:dim(MSe)[1]){
  darkmode()
  print(image(Brain2, mz = MSe$mz[i], smooth.image = "gaussian", 
              plusminus = 0.003, colorscale=magma, 
              contrast.enhance="suppression", normalize.image = "linear"))
  legend("topleft", legend= paste("correlation = ", round(MSe$correlation[i], 2)))
}
dev.off()

#(5.2) PC38:6
mycoloc <- as.data.frame(coloc_844)
MSe <- mycoloc[mycoloc$correlation >= 0.90,]
pdf(file = file.path("Result/colocalization_844.pdf"), onefile = TRUE)
for(i in 1:dim(MSe)[1]){
  darkmode()
  print(image(Brain2, mz = MSe$mz[i], smooth.image = "gaussian", 
              plusminus = 0.003, colorscale=magma, 
              contrast.enhance="suppression", normalize.image = "linear"))
  legend("topleft", legend= paste("correlation = ", round(MSe$correlation[i], 2)))
}
dev.off()

#(5.3) PC38:6
mycoloc <- as.data.frame(coloc_872)
MSe <- mycoloc[mycoloc$correlation >= 0.90,]
pdf(file = file.path("Result/colocalization_872.pdf"), onefile = TRUE)
for(i in 1:dim(MSe)[1]){
  darkmode()
  print(image(Brain2, mz = MSe$mz[i], smooth.image = "gaussian", 
              plusminus = 0.003, colorscale=magma, 
              contrast.enhance="suppression", normalize.image = "linear"))
  legend("topleft", legend= paste("correlation = ", round(MSe$correlation[i], 2)))
}
dev.off()

Session information

sessioninfo::session_info()
## ─ Session info ───────────────────────────────────────────────────────────────
##  setting  value                       
##  version  R version 4.0.2 (2020-06-22)
##  os       macOS Catalina 10.15.7      
##  system   x86_64, darwin17.0          
##  ui       X11                         
##  language (EN)                        
##  collate  en_US.UTF-8                 
##  ctype    en_US.UTF-8                 
##  tz       Asia/Jerusalem              
##  date     2020-10-19                  
## 
## ─ Packages ───────────────────────────────────────────────────────────────────
##  package      * version  date       lib source        
##  abind          1.4-5    2016-07-21 [1] CRAN (R 4.0.0)
##  assertthat     0.2.1    2019-03-21 [1] CRAN (R 4.0.0)
##  biglm          0.9-2    2020-05-26 [1] CRAN (R 4.0.0)
##  Biobase        2.48.0   2020-04-27 [1] Bioconductor  
##  BiocGenerics * 0.34.0   2020-04-27 [1] Bioconductor  
##  BiocParallel * 1.22.0   2020-04-27 [1] Bioconductor  
##  bitops         1.0-6    2013-08-17 [1] CRAN (R 4.0.0)
##  Cardinal     * 2.6.0    2020-04-27 [1] Bioconductor  
##  cli            2.1.0    2020-10-12 [1] CRAN (R 4.0.2)
##  colorspace     1.4-1    2019-03-18 [1] CRAN (R 4.0.0)
##  crayon         1.3.4    2017-09-16 [1] CRAN (R 4.0.0)
##  crosstalk      1.1.0.1  2020-03-13 [1] CRAN (R 4.0.0)
##  data.table     1.13.0   2020-07-24 [1] CRAN (R 4.0.2)
##  DBI            1.1.0    2019-12-15 [1] CRAN (R 4.0.0)
##  digest         0.6.25   2020-02-23 [1] CRAN (R 4.0.0)
##  dplyr          1.0.2    2020-08-18 [1] CRAN (R 4.0.2)
##  EBImage      * 4.30.0   2020-04-27 [1] Bioconductor  
##  ellipsis       0.3.1    2020-05-15 [1] CRAN (R 4.0.0)
##  evaluate       0.14     2019-05-28 [1] CRAN (R 4.0.0)
##  fansi          0.4.1    2020-01-08 [1] CRAN (R 4.0.0)
##  farver         2.0.3    2020-01-16 [1] CRAN (R 4.0.0)
##  fftwtools      0.9-9    2020-10-03 [1] CRAN (R 4.0.2)
##  generics       0.0.2    2018-11-29 [1] CRAN (R 4.0.0)
##  ggplot2      * 3.3.2    2020-06-19 [1] CRAN (R 4.0.0)
##  glue           1.4.2    2020-08-27 [1] CRAN (R 4.0.2)
##  gtable         0.3.0    2019-03-25 [1] CRAN (R 4.0.0)
##  htmltools      0.5.0    2020-06-16 [1] CRAN (R 4.0.0)
##  htmlwidgets    1.5.2    2020-10-03 [1] CRAN (R 4.0.2)
##  httr           1.4.2    2020-07-20 [1] CRAN (R 4.0.2)
##  irlba          2.3.3    2019-02-05 [1] CRAN (R 4.0.0)
##  jpeg           0.1-8.1  2019-10-24 [1] CRAN (R 4.0.0)
##  jsonlite       1.7.1    2020-09-07 [1] CRAN (R 4.0.2)
##  knitr          1.30     2020-09-22 [1] CRAN (R 4.0.2)
##  labeling       0.3      2014-08-23 [1] CRAN (R 4.0.0)
##  lattice        0.20-41  2020-04-02 [1] CRAN (R 4.0.2)
##  lazyeval       0.2.2    2019-03-15 [1] CRAN (R 4.0.0)
##  lifecycle      0.2.0    2020-03-06 [1] CRAN (R 4.0.0)
##  locfit         1.5-9.4  2020-03-25 [1] CRAN (R 4.0.0)
##  magrittr       1.5      2014-11-22 [1] CRAN (R 4.0.0)
##  MASS           7.3-53   2020-09-09 [1] CRAN (R 4.0.2)
##  Matrix         1.2-18   2019-11-27 [1] CRAN (R 4.0.2)
##  matter         1.14.0   2020-04-27 [1] Bioconductor  
##  mclust         5.4.6    2020-04-11 [1] CRAN (R 4.0.0)
##  munsell        0.5.0    2018-06-12 [1] CRAN (R 4.0.0)
##  nlme           3.1-149  2020-08-23 [1] CRAN (R 4.0.2)
##  pillar         1.4.6    2020-07-10 [1] CRAN (R 4.0.2)
##  pkgconfig      2.0.3    2019-09-22 [1] CRAN (R 4.0.0)
##  plotly       * 4.9.2.1  2020-04-04 [1] CRAN (R 4.0.2)
##  png            0.1-7    2013-12-03 [1] CRAN (R 4.0.0)
##  ProtGenerics * 1.20.0   2020-04-27 [1] Bioconductor  
##  purrr          0.3.4    2020-04-17 [1] CRAN (R 4.0.0)
##  R6             2.4.1    2019-11-12 [1] CRAN (R 4.0.0)
##  Rcpp           1.0.5    2020-07-06 [1] CRAN (R 4.0.2)
##  RCurl          1.98-1.2 2020-04-18 [1] CRAN (R 4.0.0)
##  rlang          0.4.8    2020-10-08 [1] CRAN (R 4.0.2)
##  rmarkdown      2.4      2020-09-30 [1] CRAN (R 4.0.2)
##  S4Vectors    * 0.26.1   2020-05-16 [1] Bioconductor  
##  scales         1.1.1    2020-05-11 [1] CRAN (R 4.0.0)
##  sessioninfo    1.1.1    2018-11-05 [1] CRAN (R 4.0.2)
##  signal         0.7-6    2015-07-30 [1] CRAN (R 4.0.0)
##  sp             1.4-4    2020-10-07 [1] CRAN (R 4.0.2)
##  stringi        1.5.3    2020-09-09 [1] CRAN (R 4.0.2)
##  stringr        1.4.0    2019-02-10 [1] CRAN (R 4.0.0)
##  tibble         3.0.3    2020-07-10 [1] CRAN (R 4.0.2)
##  tidyr          1.1.2    2020-08-27 [1] CRAN (R 4.0.2)
##  tidyselect     1.1.0    2020-05-11 [1] CRAN (R 4.0.0)
##  tiff           0.1-5    2013-09-04 [1] CRAN (R 4.0.0)
##  vctrs          0.3.4    2020-08-29 [1] CRAN (R 4.0.2)
##  viridisLite    0.3.0    2018-02-01 [1] CRAN (R 4.0.0)
##  withr          2.3.0    2020-09-22 [1] CRAN (R 4.0.2)
##  xfun           0.18     2020-09-29 [1] CRAN (R 4.0.2)
##  yaml           2.2.1    2020-02-01 [1] CRAN (R 4.0.0)
## 
## [1] /Library/Frameworks/R.framework/Versions/4.0/Resources/library